Single-cell transcriptomics uncovers EGFR signaling-mediated gastric progenitor cell differentiation in stomach homeostasis

Defects in gastric progenitor cell differentiation are associated with various gastric disorders, including atrophic gastritis, intestinal metaplasia, and gastric cancer. However, the mechanisms underlying the multilineage differentiation of gastric progenitor cells during healthy homeostasis remain poorly understood. Here, using a single-cell RNA sequencing method, Quartz-Seq2, we analyzed the gene expression dynamics of progenitor cell differentiation toward pit cell, neck cell, and parietal cell lineages in healthy adult mouse corpus tissues. Enrichment analysis of pseudotime-dependent genes and a gastric organoid assay revealed that EGFR-ERK signaling promotes pit cell differentiation, whereas NF-κB signaling maintains gastric progenitor cells in an undifferentiated state. In addition, pharmacological inhibition of EGFR in vivo resulted in a decreased number of pit cells. Although activation of EGFR signaling in gastric progenitor cells has been suggested as one of the major inducers of gastric cancers, our findings unexpectedly identified that EGFR signaling exerts a differentiation-promoting function, not a mitogenic function, in normal gastric homeostasis.

Supplementary Fig. 1 Single cell sorting of gastric cells. Related to Fig. 1 a Single cell soring of dataset 1. After the removal of small events (left, R1), dead cells were removed based on PI staining (middle, R7), and the cells in the singlet gate (right, R1&R7&R9) were sorted to 384 well plates.
b Single cell sorting of dataset 2. After the removal of dead cells by TO-PRO-3 staining (left, R9), the cells with nuclei were selected based on Hoechst 33342 (middle, R7), and cells in the singlet gate (right, R9&R7&R5) were sorted to 384 well plates.
c Phase contrast images of dissociated gastric cells used for single cell sorting in dataset 2. Before cell sorting (left), the cell sample contains debris (white arrowheads). The singlet gate R9&R7&R5 in (b) contains neither debris nor cell aggregates (middle). The aggregate gate R9&R7&R8 in (b) contains cell doublets, triplets (black arrowheads, right), and some larger cells. Scale bar, 50 μm. These images are representative of three images from one single cell sorting experiment. 10,000 1 3 4 9 3 6 5 1 1 0 2 7 1 2 1 0 8 1 1 4 1 5 1 6 1 7 1 8 1 9 Mitochondrial ratio per cell (%) Mitochondrial ratio per cell (%) Gene count per cell  Supplementary Fig. 2 Quality control of single-cell RNA-seq data. Related to Fig. 1 a Scatter plot showing gene count and mitochondrial ratio of all sequenced cells in dataset 1. The cells are colored by plate ID. The cells with less than 2,500 gene counts were removed as low-quality cells.
b Violin plots showing gene count and mitochondrial ratio of the cells that passed quality control in each cluster. Two biological replicates prepared by different operators were shown in red and green.
c Scatter plot showing gene count and mitochondrial ratio of all sequenced cells in dataset 2. The cells are colored by plate ID. The cells with less than 2,500 gene counts were removed as low-quality cells.
d Violin plots showing gene count and mitochondrial ratio of the cells that passed quality control in each cluster. Two biological replicates prepared by different operators were shown in red and green. e UMAP of gastric cells in dataset 1 (left, our study), dataset 2 (middle, our study), and the previously published dataset (right, GSE157694). In our study, the cells were colored by plate ID. In the previously published data, the cells were colored by plates. Note that, in our study, the cells from six plates (dataset 1, two biological replicates) or twelve plates (dataset 2, two biological replicates) were clustered together, suggesting that our data is not affected by technical noise.   Supplementary Fig. 3. Single-cell analysis of corpus and antrum gastric units (dataset 1).
a UMAP visualization of the cells isolated from the corpus and antrum gastric units (dataset 1). Cells are colored according to the clustering results.
b Violin plots of dataset 1 showing the expression of known marker genes in each cluster.
c Western blotting showing the expression of BASP1 (top), GKN2 (middle), and α-Tubulin (bottom) in corpus and antrum glands isolated from three different mice. Data are presented as mean relative expression ± SD. Significance was calculated by two-tailed Student's t-tests (**p = 0.0042). Source data are provided as a Source Data file.
d Immunofluorescence staining of adult mouse corpus and antrum tissue with KRT7 (green), KRT20 (green), MKI67 (red), and DAPI (blue). Scale bar, 100 μm. The images of KRT20, GKN2, and AQP3 staining are representative of three independent experiments. The images of KRT7 staining are representative of two independent experiments.   a Immunofluorescence staining of the corpus tissues with cell type-specific markers. High magnification images of the dotted squares are shown on the right side. All images are representative of two independent experiments.
b Violin plots of dataset 2 showing the expression of previously identified gastric stem cell markers, isthmus progenitor cell-enriched markers identified in this study, and the markers that distinguish between corpus and antrum pit cells.

Muc5ac Pgc Muc5ac_Pgc
Pgc 4  b UMAP visualization of isthmus progenitor cells identified in pit cell or neck cell lineage by FateID analysis. Note that isthmus progenitor cells used for pit cell (yellow) or neck cell (blue) lineage by FateID were merged with isthmus progenitor cells with high pit cell feature or high neck cell feature identified in Fig. 1f, supporting the reliability of FateID analysis.  Supplementary Fig. 9. Characterization of the pit cell differentiation process in dataset 1. Related to Fig.2. a Pseudotemporal ordering of major corpus epithelial cells identified in dataset 1. Principal curves are shown for mature pit cells (t8, heavy dotted line), neck cells (t3, solid line), and parietal cells (t99, fine dotted line).
b Self-organizing map of binarized pseudotemporal expression profiles along the pit cell differentiation trajectory. The xaxis indicates the cells involved in the pit cell lineage that are colored according to cell types and the y-axis indicates the nodes.
c Hierarchical clustering of co-expression nodes in (b). The nodes in group A and B show upregulation and downregulation patterns of the genes in the pit cell differentiation trajectory, respectively. The x-axis indicates the cells involved in the pit cell lineage that are colored according to cell types and the y-axis indicates the nodes. The x-axis in (c) is not same as that of (b) due to the new hierarchical clustering of the cells. The colors for each cell type are shown on the right bottom corner.   a Self-organizing map of binarized pseudotemporal expression profiles along the neck cell differentiation trajectory. The x-axis indicates the cells involved in the neck cell lineage that are colored according to cell types and the y-axis indicates the nodes.

GO (Cellular Component)
b Hierarchical clustering of co-expression nodes in (b). The nodes in group C and D show upregulation and downregulation patterns of the genes in the neck cell differentiation trajectory, respectively. The x-axis indicates the cells involved in the neck cell lineage that are colored according to cell types and the y-axis indicates the nodes. The x-axis in (b) is not same as that of (a) because of the new hierarchical clustering of the cells. The colors for each cell type are shown on the right bottom corner.
c Characteristics of group C and D identified in (b). Marker genes in each group are shown. Average pseudotemporal expression profile of the representative node is shown in black line. The x-axis indicates pseudotime from isthmus progenitor cells to neck cells and the y-axis indicates expression level. Colors in each dot represent cell types, which corresponds to colors used in Fig. 2a Supplementary Fig. 11. Characterization of the neck cell differentiation process in dataset 1. Related to Fig.2. a Self-organizing map of binarized pseudotemporal expression profiles along the neck cell differentiation trajectory. The x-axis indicates the cells involved in the neck cell lineage that are colored according to cell types and the y-axis indicates the nodes.
b Hierarchical clustering of co-expression nodes in (b). The nodes in group C and D show upregulation and downregulation patterns of the genes in the neck cell differentiation trajectory, respectively. The x-axis indicates the cells involved in the neck cell lineage that are colored according to cell types and the y-axis indicates the nodes. The x-axis in (b) is not same as that of (a) because of the new hierarchical clustering of the cells. The colors for each cell type are shown on the right bottom corner.
c Characteristics of group C and D identified in (b). Previously reported marker genes in each group are shown. Average pseudotemporal expression profile of the representative node is shown in gray line. The x-axis indicates pseudotime from isthmus progenitor cells to neck cells and the y-axis indicates expression level. Colors in each dot represent cell types, which corresponds to colors used in Supplementary Fig. 7a  a UMAP visualization of gastric cells other than the major corpus epithelial identified in datasets 1 and 2.
b Violin plots showing mitochondrial ratio per cell in each cluster. The cells derived from datasets 1 and 2 are in red and green, respectively.
c Violin plots showing the expression of known marker genes in each cluster.
d Immunofluorescence staining of adult mouse corpus and antrum tissues with ACTA2 (green), PDGFRα (red), PECAM1 (red), CD4 (red), and DAPI (blue). High magnification images of the dotted squares are shown on the right side. As for CD4 staining, phase contrast pictures merged with CD4 and DAPI are shown on the right side. Yellow lines indicate each gastric gland. The images are representative of three independent experiments (ACTA2/PDGFRa) or two independent experiments (ACTA2/PECAM1 and CD4).    Fig. 15. CellChat-based Analysis of outgoing and incoming signaling in each cell type. Related to Fig. 3 The x-axis indicates cell types and the y-axis indicates signaling pathways. For each signaling, the outgoing signal intensity is shown in red (upper row) and the incoming signal intensity is shown in blue (lower row).  a Immunofluorescence staining of MKI67 (green), IGF1R (red), and DAPI (blue) in adult mouse corpus tissues. High magnification images of the dotted squares are shown on the right side. The images are representative of two independent experiments.
b Relative contribution of BMP, WNT, SHH, and IGF ligand-receptor pairs calculated by CellChat.
c Violin plots showing the expression of IGF, BMP, SHH, and WNT ligands and receptors that were identified as major contributors by CellChat.   Supplementary Fig. 19. Expression pattern of ligands, receptors, and downstream regulators for IGF (a), NOTCH (b), BMP (c), GDF (d), and Interleukin (e) signaling pathways across all gastric cell types. a qPCR analysis of gastric epithelial cell marker expression in erlotinib-treated corpus organoids in Fig. 4f. Data are presented as mean fold changes ± SD (n = 3 biologically independent samples for each culture condition). **p = 0.0005 for Stmn1; **p = 0.0093 for Tff2; *p = 0.0131 for Troy; *p = 0.0179 for Gif.
b Immunofluorescence analysis of pERK (red) and DAPI (blue) in corpus organoids cultured with or without TGFα and erlotinib. Scale bar, 200 µm.
c Quantification of pERK fluorescence intensity in corpus organoids in (b). Data are presented as mean fold changes ± SD (n = 3 biologically independent samples). Each data point represents the mean value of at least ten organoids. TGFa vs Control p < 0.0001, TGFa + erlotinib vs Control n.s., TGFa + erlotinib vs TGFa p < 0.0001.
d Immunofluorescence staining of GKN2 (green), MKI67 (red), SOX9 (green), PGC (green), and DAPI (blue) in corpus organoids treated with or without erlotinib (0.5 µM). The adult mouse corpus tissues stained with PGC (green) is shown in the upper right corner as a positive control. Scale bar, 200 µm.
e Quantification of GKN2 fluorescence intensity and the percentage of MKI67 + and SOX9 + cells in the corpus organoids in (d). Data are presented as mean fold changes ± SD (n = 3 or 4 biologically independent samples). **p = 0.0049 for GKN2 fluorescence intensity; **p = 0.0014 for SOX9 + nuclei.
f The control corpus organoid and the corpus tissue stained with anti-AQP3 antibody (green), phalloidin (red), and DAPI (blue). Images were captured by laser confocal microscopy and representative of three independent experiments.
g Comparison of growth factors and chemicals used in the previous study (Barker et al. 2010) and those used in our study.
h qPCR analysis of gastric epithelial cell marker expression in the corpus organoids cultured with rWnt3a (100 ng/mL) or Wnt3a CM. The response to EGF (50 ng/mL) was analyzed. Data are presented as mean fold changes ± SD (n = 3 biologically independent samples for each culture condition  a qPCR analysis of gastric epithelial cell markers in the corpus gastric organoids cultured with or without 100 ng/mL TNFSF12 in Fig. 5a. Data are presented as mean fold changes ± SD (n = 3 biologically independent samples for each culture condition). **p = 0.0009 for Stmn1; *p = 0.0390 for Tnfrsf12a; *p = 0.0164 for Tff2; **p = 0.0009 for Troy.
c Quantification of the fluorescence intensity of NF-κB1 in the corpus organoids in (b). Data are presented as mean fold changes ± SD (n = 3 biologically independent samples). Each data point represents the mean value of at least ten organoids. TGFa vs Control p = 0.0485, TNFSF12 vs Control p = 0.0006, TNFSF12 vs TGFa p < 0.0001.
d qPCR analysis of gastric epithelial cell markers in corpus gastric organoids cultured with or without 100 ng/mL TNFα in Fig. 5e. Data are presented as mean fold changes ± SD (n = 3 biologically independent samples for each culture condition). *p = 0.0139 for Stmn1; *p = 0.0393 for Tff2.
f Quantification of the fluorescence intensity of NF-κB1 in the gastric organoids in (e). Data are presented as mean fold changes ± SD (n = 3 biologically independent samples). Each data point represents the mean value of ten organoids. **p = 0.0066.
h Quantification of GKN2 fluorescence intensity and the percentage of MKI67 + and SOX9 + cells of the corpus gastric organoids in (g). Data are presented as mean fold changes ± SD (n = 3 biologically independent samples). Each data point represents the mean value of at least ten organoids. *p = 0.0220 for GKN2 fluorescence intensity, *p = 0.0126 for MKI67 + nuclei.
Statistical information: significance was calculated by two-tailed Student's t-tests for samples with equal variances or two-sided Welch's t-tests for samples with unequal variances in (a), (d), (f), and (h)(*p < 0.05; ** p < 0.01); significance was calculated by one-way ANOVA followed by Tukey post hoc test at the 0.05 significance level in (c). Source data are provided as a Source Data file. a qPCR analysis of gastric epithelial cell markers in corpus gastric organoids cultured with or without a NF-κB inhibitor, QNZ (16 nM) in Fig. 6a. Data are presented as mean fold changes ± SD (n = 3 biologically independent samples). ** p = 0.0036 for Tff2.
b qPCR analysis of gastric epithelial cell marker expression in corpus gastric organoids cultured with or without a NF-κB inhibitor NaSal (1.25 mM) in Fig. 6c. Data are presented as mean fold changes ± SD (n = 3 biologically independent samples).
e Immunofluorescence staining of PGC (green) and DAPI (blue) in the corpus gastric organoids cultured with or without 16 nM QNZ or 1.25 mM NaSal. Scale bar, 200 µm. The images are representative of two biologically independent samples.
g qPCR analysis of gastric epithelial cell marker expression in the corpus organoids in (f). Data are presented as mean fold changes ± SD (n = 4 biologically independent samples Statistical information: significance was calculated by two-tailed Student's t-tests for samples with equal variances or two-sided Welch's t-tests for samples with unequal variances in (a), (b), (d). * p < 0.05; ** p < 0.01. In (g) and (i), significance was calculated by one-way ANOVA followed by Tukey post hoc test (significance level, p < 0.05). Source data are provided as a Source Data file.  Supplementary Fig. 25. Activation of ERK pathway is required for pit cell differentiation. Related to Fig. 7 a qPCR analysis of gastric epithelial cell markers in corpus organoids cultured with or without TGFα (12.5 ng/mL) and an ERK1/2 inhibitor, SCH772984 (1 μM). Data are presented as mean fold changes ± SD (n = 3 biologically independent samples). TGFa vs Control p = 0.0197, TGFa + SCH vs Control p = 0.0449, TGFa + SCH vs TGFa p = 0.0010 for Stmn1; p = 0.0437, n.s., p = 0.0196 for Tff2; n.s., n.s., p = 0.0072 for Troy; n.s., n.s., n.s. for Lgr5.
c Corpus organoids treated with or without a Stat3 inhibitor Stattic. The organoids were cultured in the presence of TGFα (12.5 ng/mL) with or without Stattic (50 μM) for 6 days. Scale bar, 500 µm.
d qPCR analysis of gastric epithelial cell marker expression in the corpus organoids in (c). Data are presented as mean fold changes ± SD (n = 3 biologically independent samples for each culture condition). Note that Stattic did not affect the expression of pit cell markers and Mki67. e Immunofluorescence staining of pERK (red), MKI67 (green), and DAPI (blue) in the mouse antrum tissue. Inset shows the high magnification image of boxed area. The image are a representative of three independent experiments.
Statistical information: significance was calculated by one-way ANOVA followed by Tukey post hoc test in (a), (b), and (d) (significance level, p < 0.05).